{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Lecture 3 : Numbers and symbolic expressions\n",
    "**References:**\n",
    "* [Programming in Python and Sage](https://doc.sagemath.org/html/en/thematic_tutorials/tutorial-programming-python.html)\n",
    "* [Introduction to Symbolic Computation](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/index.html), Lectures [4](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec04.html), [5](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec05.html), [6](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec06.html), [7](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec07.html),\n",
    "[14](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec14.html)\n",
    "[28](http://homepages.math.uic.edu/~jan/mcs320/mcs320notes/lec28.html)\n",
    "* [Basic properties of the integers in Sage](https://doc.sagemath.org/html/en/thematic_tutorials/group_theory.html#basic-properties-of-the-integers)\n",
    "\n",
    "**Summary:**<br>\n",
    "In todays lecture, we'll start with a reminder of **conditional execution** (``if``) and **loop** (``for``) in Python, and the datatype of **dictionaries**. Then we discuss various types of numbers in SageMath, such as **integers**, **rationals**, **algebraic numbers** and more. In the last part, we discuss **symbolic expressions** and how to use them to solve equations."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Python warmup : conditional execution (if/else)\n",
    "Sometimes we want to perform some operations in our code (like computing the inverse of a matrix $A$) only if certain conditions are satisfied (the matrix $A$ is invertible), and maybe do something else otherwise (print a warning). This can be done in Python using the if/else construction:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = matrix([[1,2],[3,5]]); A"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "if A.is_invertible():\n",
    "    print('The inverse of A is:')\n",
    "    print(A.inverse())\n",
    "else:\n",
    "    print('A is not invertible!')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here the general pattern is:\n",
    "```\n",
    "if CONDITION:\n",
    "    do_something()\n",
    "    do_more_of_something()\n",
    "else:\n",
    "    do_something_else()\n",
    "    maybe_more_of_something_else()\n",
    "```\n",
    "Note that the indentation (using a ``Tab`` = 4 spaces) is important!\n",
    "\n",
    "Here ``CONDITION`` is an expressiong evaluating to a ``bool``, which is either ``True`` or ``False``. Basic ingredients for such conditions:\n",
    "* ``A == B``, ``A != B`` for checking equality (or not-equality)\n",
    "* ``A <= B``, ``A > B`` etc for comparisons of two values\n",
    "* Logical operations: ``Cond1 and Cond2``, ``Cond1 or Cond2``, ``not Cond1``\n",
    "\n",
    "It is possible to omit the ``else`` part, in which case nothing is done if ``CONDITION`` is not satisfied."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Given real parameters $p,q$ of the quadratic equation\n",
    "$$\n",
    "x^2 + p x + q = 0\n",
    "$$\n",
    "write some code which prints the set of solutions $x \\in \\mathbb{R}$ (which you have to compute yourself, since solving equations only comes later in the lecture).\n",
    "\n",
    "*Hint:* It's never too late in life to look at the [formula for the solutions of a general quadratic equation](https://en.wikipedia.org/wiki/Quadratic_formula)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution:** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "if p^2/4-q > 0:\n",
    "    print('The equation has two solutions:')\n",
    "    print(-p/2 + sqrt(p^2/4-q))\n",
    "    print(-p/2 - sqrt(p^2/4-q))\n",
    "if p^2/4-q == 0:\n",
    "    print('The equation has one solution:')\n",
    "    print(-p/2)\n",
    "if p^2/4-q < 0:\n",
    "    print('The equation has no solutions.')\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that we can have an ``if``-block inside another ``if``- or ``else``-block, using double indentation (and so on):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "A = -6\n",
    "\n",
    "if A.is_prime():\n",
    "    if A == 2:\n",
    "        print('A is an even prime')\n",
    "    else:\n",
    "        print('A is an odd prime')\n",
    "else:\n",
    "    if A < 0 and (-A).is_prime():\n",
    "        print('A is negative, but -A is a prime.')\n",
    "    else:\n",
    "        print('A is not a prime, since it factors as:')\n",
    "        print(factor(A))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Python warmup: for-loops\n",
    "Sometimes we want to perform an operation for all $n$ in a range of values. This can be accomplished by a ``for``-loop. As an example, below we compute the sum of the numbers $n=0,1, \\ldots, 100$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "count = 0\n",
    "for n in range(101):\n",
    "    count = count + n\n",
    "print(count)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The general pattern is:\n",
    "```\n",
    "for n in ITERABLE:\n",
    "    do_something()\n",
    "```\n",
    "Above, we have used ``range(N)`` as the ``ITERABLE``, which outputs the numbers $0,1, \\ldots, N-1$, which is why we had to take ``range(101)`` to sum from $0$ to $100$. "
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Count for how many of the numbers $n=0,1, \\ldots, 2022$ the value $\\sin(n)$ is negative."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)<br>\n",
    "<!---\n",
    "```\n",
    "count = 0\n",
    "for n in range(2023):\n",
    "    if sin(n)<0:\n",
    "        count += 1\n",
    "count\n",
    "> 1011\n",
    "```\n",
    "<img src=\"https://www.meme-arsenal.com/memes/6b90c98fda17445646a21305bc75bfc7.jpg\" width=300/>\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Similarly, we have the variations:\n",
    "* ``range(N0, N)`` for the numbers $N_0, N_0+1, \\ldots, N-1$,\n",
    "* ``range(N0, N, d)`` for the numbers $N_0, N_0+d, N_0+2d, \\ldots, N_0+k\\cdot d$, where $k$ is maximal such that $N_0+k\\cdot d<N$\n",
    "\n",
    "Some more interesting examples:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for a in [2,4,'cat']:\n",
    "    print(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we are not just using Python, but SageMath, we also have some cool iterables, like the symmetric group $S_n$ given by ``SymmetricGroup(n)``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "for a in SymmetricGroup(3):\n",
    "    print('The element {} of S_3 has order {}'.format(a, a.order()))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For more information about iterables, see [this](https://doc.sagemath.org/html/en/thematic_tutorials/tutorial-comprehensions.html) tutorial."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Python warmup: dictionaries\n",
    "Dictionaries are a useful data structure in Python, which encodes a map from a finite set of objects to some other objects:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# The dictionary d sending 1 -> 2, 2 -> 2, 3 -> 'cat'\n",
    "d = {1:2, 2:2, 3: 'cat'}"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d[1]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can easily add new assignments to the dictionary (or change old ones):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d[4] = 3\n",
    "d[3] = 'dog'\n",
    "d"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "There is also a connection to the ``for``-loops we saw above: on the one hand, dictionaries are themselves ``Iterables``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "for a in d:\n",
    "    print(a,d[a])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "On the other hand, similar to the constructions of lists like ``[k^2 for k in range(10)]`` that we already saw, we can also create dictionaries in this way:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "s = {i : i+1 for i in range(10)}\n",
    "s"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "* Create dictionaries ``f,g`` implementing the maps\n",
    "$$\n",
    "f : \\{1,2,3\\} \\to \\{5,6\\}, 1 \\mapsto 5, 2 \\mapsto 6, 3 \\mapsto 5,\\\\\n",
    "g : \\{4,5,6\\} \\to \\{7,8,9\\}, 4 \\mapsto 9, 5 \\mapsto 8, 6 \\mapsto 9.\n",
    "$$"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "* Create the dictionary ``h`` implementing the composition $g \\circ f$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "f = {1:5, 2:6, 3:5}\n",
    "g = {4:9, 5:8, 6:9}\n",
    "h = {a: g[f[a]] for a in f}\n",
    "h\n",
    "> {1: 8, 2: 9, 3: 8}\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Numbers\n",
    "Last time, we already saw some interesting rings implemented in SageMath:\n",
    "* ``ZZ`` the integers\n",
    "* ``QQ`` the rational numbers\n",
    "* ``RR`` the real numbers (floating point operations!)\n",
    "* ``CC`` the complex numbers (floating point operations!)\n",
    "* ``GF(q)`` the finite field $\\mathbb{F}_q$ with $q$ elements\n",
    "\n",
    "In this section, we will learn more about elements of these rings (commonly known as \"numbers\") and what we can do with them.\n",
    "\n",
    "### Integers\n",
    "When we type an expression like ``42`` in SageMath, this will produce an Integer:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(42)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Integer division and remainder\n",
    "Apart from the usual division ``a/b`` of two Integers (which mostly produces a Rational), we can also use ``a//b`` to compute the largest integer which is at most $a/b$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "5/3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "5//3"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "6//3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The remainder (or modulus) of this division can be computed with ``a%b``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "5 % 3"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Check whether the number $17^{17}-1$ is divisible by $3$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "(17^17-1) % 3\n",
    "> 1\n",
    "```\n",
    "Thus the number is not divisible by $3$, it has remainder $1$ when dividing it by $3$.\n",
    "--->"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Compute the last $3$ digits of the integer\n",
    "$$\n",
    "a = (((2^2)^2)^\\ldots)^2\n",
    "$$\n",
    "where the number $2$ appears \n",
    "* $5$ times, or\n",
    "* $500$ times \n",
    "\n",
    "in the formula for $a$.<br><br>\n",
    "*Hint 1:* Since the number of digits of an integer roughly doubles when taking its square, we can estimate that $a$ has more decimal digits than the [number of atoms in the universe](https://en.wikipedia.org/wiki/Eddington_number). Keep in mind that the code you use for solving the exercise must run on a computer which is part of the universe!<br>\n",
    "*Hint 2:* If you have ignored Hint 1, you probably at some point want to click on ``Kernel -> Interrupt`` in the menu to stop the computation, and maybe also on ``Kernel -> Restart`` to clear your memory."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)<br>\n",
    "<!---\n",
    "To get the last three digits of a number ``a`` we can compute ``a % 1000``. On the other hand, we know from elementary number theory that the last three digits of the square of ``a`` only depend on the last three digits of ``a``:\n",
    "```\n",
    "a = 5276582\n",
    "a^2\n",
    "> 27842317602724\n",
    "(a % 1000)^2\n",
    "> 338724\n",
    "```\n",
    "So we can just start with ``a = 2`` and then 499 times\n",
    "* take the square, and\n",
    "* reduce modulo 1000\n",
    "\n",
    "for which we can use a ``for``-loop:\n",
    "```\n",
    "a = 2\n",
    "for i in range(499):\n",
    "    a = a^2 % 1000\n",
    "a\n",
    "> 56\n",
    "```\n",
    "Thus the last three digits are $056$.\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Factorization\n",
    "We can also compute the prime factorization of an integer:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "u = factor(99); u"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(u)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To actually extract the information from this factorization, to use it in further computations, the following conversions are useful:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "list(u)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "dict(u)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Find the largest $m \\geq 1$ such that the number $100!$ is divisible by $2^m$."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "u = factor(factorial(100))\n",
    "u\n",
    "> 2^97 * 3^48 * 5^24 * 7^16 * 11^9 * 13^7 * 17^5 * 19^5 * 23^4 * 29^3 * 31^3 * 37^2 * 41^2 * 43^2 * 47^2 * 53 * 59 * 61 * 67 * 71 * 73 * 79 * 83 * 89 * 97\n",
    "dict(u)[2]\n",
    "> 97\n",
    "```\n",
    "Thus $m$ is equal to $97$.\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Conversely, we can look up the primes in the interval $[a,b-1]$ using ``primes(a,b)``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "primes(4,19)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Oops, this function returns a **generator**, i.e. an iterable that can be used in a ``for`` loop, without having to first compute all the numbers it will put out:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for p in primes(1,10^100):\n",
    "    print(p)\n",
    "    sleep(1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To just look at the numbers, you can convert this generator object into a list:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "list(primes(4,19))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Python int vs. SageMath Integer\n",
    "Remember the exercise from the first lecture to create a multiplication table of the integers from $1$ to $10$:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "matrix([[a*b for a in range(1,11)] for b in range(1,11)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's see what happens instead if we wanted to create a *division* table:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "matrix([[a/b for a in range(1,11)] for b in range(1,11)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "By golly, what happened? It turns out that the integers that ``range`` spits out are actually of type ``int`` in Python, and not of type ``Integer``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for i in range(10):\n",
    "    print(type(i))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Therefore, the division above is a division of ``int`` values, which returns a ``float``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = int(3)/int(7); a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To repair this, you can convert either numerator or denominator of the fraction into an ``Integer``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = ZZ(int(3))/int(7); b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Thus a nice and mathematically accurate division table can be obtained as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "matrix([[ZZ(a)/b for a in range(1,11)] for b in range(1,11)])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This type of conversion works more generally:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "c = 6/2\n",
    "(c,type(c))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "d = ZZ(c)\n",
    "(d,type(d))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "An integer $a$ is a perfect power if there exist integers $b,e$ with $e>1$ such that $a=b^e$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "[a.is_perfect_power() for a in [3, 9, 6]]"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the integers $1, \\ldots, 20$ print line by line whether they are perfect powers or not, e.g. the output should start like this:\n",
    "```\n",
    "1 True\n",
    "2 False\n",
    "3 False\n",
    "4 True\n",
    "5 False\n",
    "...\n",
    "```"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "for a in range(1,21):\n",
    "    print(a, ZZ(a).is_perfect_power())\n",
    "> \n",
    "1 True\n",
    "2 False\n",
    "3 False\n",
    "4 True\n",
    "5 False\n",
    "6 False\n",
    "7 False\n",
    "8 True\n",
    "9 True\n",
    "10 False\n",
    "11 False\n",
    "12 False\n",
    "13 False\n",
    "14 False\n",
    "15 False\n",
    "16 True\n",
    "17 False\n",
    "18 False\n",
    "19 False\n",
    "20 False\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Rational numbers\n",
    "If you want to do an exact computation with real numbers (e.g. in linear algebra), it is best when you can stay inside the rational numbers, since they have an exact representation as a fraction. Many operations we saw for Integers above work analogously for Rationals:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = 24/9; a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(a.numerator(), a.denominator())"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "b = QQ(1.45); b"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "factor(b)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Some more useful functions (which also work for integers or real numbers):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = - 5/3\n",
    "floor(a)  # rounding down"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ceil(a)   # rounding up"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sign(a)   # the sign"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "abs(a)    # absolute value"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Among the numbers of the form $p/q$ for $p,q \\in \\{1, \\ldots, 20\\}$, which is closest to $\\pi$?<br>\n",
    "*Hint:* Have a variable storing the best approximation you have so far, and when going through all numbers $p/q$ replace the current optimum if one of them comes closer."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "best = 1\n",
    "diff = pi - 1\n",
    "for p in range(1,21):\n",
    "    for q in range(1,21):\n",
    "        if abs(ZZ(p)/ZZ(q)-pi) < diff:\n",
    "            best = ZZ(p)/ZZ(q)\n",
    "            diff = abs(best-pi)\n",
    "print(best)\n",
    "print(n(best))\n",
    "print(n(abs(best-pi)))\n",
    "> 19/6\n",
    "3.16666666666667\n",
    "0.0250740130768734\n",
    "```\n",
    "Note that using some clever programming, you can get the result in a single line:\n",
    "```\n",
    "min((abs(ZZ(p)/q-pi), ZZ(p)/q) for p in range(1,21) for q in range(1,21))\n",
    "> (abs(-pi + 19/6), 19/6)\n",
    "```\n",
    "This uses that tuples are ordered by first comparing the first entries and then comparing further entries when the first entries are equal.\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Real and complex numbers\n",
    "We already saw that ``RR`` and ``CC`` are the fields of real and complex numbers, implemented with floating-point arithmetic. In particular, the computations here are no longer exact:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "piRR = RR(pi); piRR"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sin(2*piRR)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "a = CC(- pi * I); a"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "exp(a)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Since we mostly talk about symbolic computations in this course, we won't go deeper into this topic for now, but here are some links where you can find more information:\n",
    "* [Computing with real numbers](https://fredrikj.net/math/ejcim2020_joh_en.pdf)\n",
    "* [Fixed and\n",
    "Arbitrary Precision Numerical Fields](https://doc.sagemath.org/pdf/en/reference/rings_numerical/rings_numerical.pdf)\n",
    "* [Sage Quickstart for Numerical Analysis](https://doc.sagemath.org/html/en/prep/Quickstarts/NumAnalysis.html)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Algebraic numbers\n",
    "One nice field in which exact computations are possible is the field $\\overline{\\mathbb{Q}}$ of **algebraic numbers**. These are complex numbers which are a zero of some polynomial with rational coefficients. The unique such polynomial of minimal degree and with leading coefficient $1$ is called the **minimal polynomial**."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "z = QQbar(3/2 + I); z"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "z.minpoly()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "z^2 - 3*z + 13/4"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can do things like taking roots of an algebraic number. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "u = sqrt(z); u"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Even though below it is displayed as a decimal number, behind the scenes it still remembers its exact representation:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "u.minpoly() # for some reason this produces an error on the SageMath installation at the University"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can do fun things, like compute an expansion of a real algebraic number as a continued fraction:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "w = sqrt(QQbar(2)); w"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "continued_fraction(w)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This says that\n",
    "$$\n",
    "\\sqrt{2} = 1 + \\frac{1}{2 + \\frac{1}{2 + \\frac{1}{2+ \\ldots} }}\n",
    "$$\n",
    "If you haven't seen it, it's a fun exercise to prove this!"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Compute the minimal polynomial and continued fraction of the golden ratio\n",
    "$$\n",
    "\\phi = \\frac{1 + \\sqrt{5}}{2}\\,.\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "phi = (1 +sqrt(5))/2\n",
    "phi.minpoly()\n",
    "> x^2 - x - 1\n",
    "continued_fraction(phi)\n",
    "> [1; 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, ...]\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Finite fields\n",
    "Given a prime number $p$, there exists the **finite field**\n",
    "$$\\mathbb{F}_p = \\mathbb{Z}/p\\mathbb{Z}$$\n",
    "of integers modulo $p$, which has precisely $p$ elements. More generally, for a prime power $q=p^r$, $(r \\geq 1)$, there exists a unique field $\\mathbb{F}_q$ with $q$ elements. It is *not* given by $\\mathbb{Z}/q \\mathbb{Z}$ though! Instead it is given as an extension\n",
    "$$\n",
    "\\mathbb{F}_q = \\mathbb{F}_p[z]/(f(z)),\n",
    "$$\n",
    "where $z$ is a variable and $f \\in \\mathbb{F}_p[z]$ is an irreducible polynomial of degree $r$.\n",
    "\n",
    "This finite field can be created in SageMath using ``FiniteField(q)`` (or ``GF(q)`` for short)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "GF(9)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "If we actually need to write elements of $\\mathbb{F}_q$, we can use the following line to specify the variable used in the definition:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "R.<z> = GF(9)\n",
    "R"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can check e.g. that $z+z+z=0$, since $z+z+z = 3 \\cdot z = 0 \\cdot z = 0$, as $3 = 0 \\in \\mathbb{F}_{3^2}$."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "print(z+z)\n",
    "print(z+z+z)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also type formulas using standard integers like ``17``. Of course by default they give elements of type ``Integer``, but SageMath is smart enough to convert them to elements of $\\mathbb{F}_p \\subseteq \\mathbb{F}_q$. We are later going to see more about how this works."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "u = z + 17; u"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "u+1"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And we can check which irreducible polynomial $f$ was used to define $\\mathbb{F}_q = \\mathbb{F}_p[z]/(f(z))$; SageMath automatically picks one for us:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "z.minimal_polynomial()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here is a cool exercise using finite fields, which hints at some quite deep mathematics:\n",
    "#### Exercise\n",
    "Given a prime power $q$, let $\\mathbb{L}_q(f)$ be the set of solutions $(x,y,z) \\in \\mathbb{F}_q^3$ of the equation\n",
    "$$\n",
    "f = x^3 y - 3x^2 z^2 - yz + 4 \\stackrel{\\text{!}}{=} 0 \\in \\mathbb{F}_q.\n",
    "$$\n",
    "* Compute the cardinality $|\\mathbb{L}_2(f)|$ by going through all points in $\\mathbb{F}_2^3$ and counting how many of them satisfy the equation.\n",
    "* Compute the cardinality $|\\mathbb{L}_q(f)|$ for a few more values of $q$ and then make a guess how it grows with $q$.<br>\n",
    "*Hint:* There is a function in SageMath giving you a list of all prime powers $q$ up to $n-1$. Can you guess its name and confirm with Tab-completion?"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "<b>Solution:</b> (uncomment to see)<br>\n",
    "<!---\n",
    "We get a list of all prime powers $q$ up to $n-1$ by ``prime_powers(n)``. To count the solutions, we can just try all values $(x,y,z) \\in \\mathbb{F}_q^3$. Finally, a convenient way to store the answer is to use a use a dictionary, sending $q$ to $|\\mathbb{L}_q(f)|$. Then we obtain something as follows:\n",
    "```\n",
    "countdict = {}\n",
    "for q in prime_powers(30):\n",
    "    count = 0\n",
    "    for x in GF(q):\n",
    "        for y in GF(q):\n",
    "            for z in GF(q):\n",
    "                if x^3*y - 3*x^2 * z^2 - x*(y*z + 4) == 0:\n",
    "                    count += 1\n",
    "    countdict[q] = count\n",
    "```\n",
    "Alternatively, we could store our results in a list ``countlist``, which we create with ``countlist = []`` and where we write ``countlist.append((q,count))`` in the last line.\n",
    "\n",
    "When we output the value of ``countdict`` we obtain:\n",
    "```\n",
    "{2: 4,\n",
    " 3: 6,\n",
    " 4: 16,\n",
    " 5: 20,\n",
    " 7: 42,\n",
    " 8: 64,\n",
    " 9: 72,\n",
    " 11: 132,\n",
    " 13: 156,\n",
    " 16: 256,\n",
    " 17: 272,\n",
    " 19: 342,\n",
    " 23: 552,\n",
    " 25: 600,\n",
    " 27: 702,\n",
    " 29: 812}\n",
    "```\n",
    "To get an idea how fast this grows, one can try to plot it (we will see the details later):\n",
    "```\n",
    "list_plot([(a,countdict[a]) for a in countdict], color='red')\n",
    "```\n",
    "<img src=\"https://johannesschmitt.gitlab.io/mat007/LangWeil1.png\" width=500/>\n",
    "\n",
    "Looking at this picture, maybe we are reminded of the normal parabola from highschool. So, if the function is roughly $q \\mapsto q^2$, then taking the square-root we should roughly get the function $q \\mapsto q$. And indeed:\n",
    "```\n",
    "{q: n(sqrt(v)) for q,v in countdict.items()}\n",
    "> {2: 2.00000000000000,\n",
    " 3: 2.44948974278318,\n",
    " 4: 4.00000000000000,\n",
    " 5: 4.47213595499958,\n",
    " 7: 6.48074069840786,\n",
    " 8: 8.00000000000000,\n",
    " 9: 8.48528137423857,\n",
    " 11: 11.4891252930761,\n",
    " 13: 12.4899959967968,\n",
    " 16: 16.0000000000000,\n",
    " 17: 16.4924225024706,\n",
    " 19: 18.4932420089069,\n",
    " 23: 23.4946802489415,\n",
    " 25: 24.4948974278318,\n",
    " 27: 26.4952825989835,\n",
    " 29: 28.4956136975500}\n",
    "```\n",
    "So we find that $|\\mathbb{L}_q(f)|\\approx q^2$. Plotting the two functions side by side, we see that they are very close:\n",
    "```\n",
    "list_plot(list(countdict.items()), color='red')+plot(lambda x: x**2,0,30)\n",
    "```\n",
    "<img src=\"https://johannesschmitt.gitlab.io/mat007/LangWeil2.png\" width=500/>\n",
    "\n",
    "Behind this phenomenon is the so-called [**Lang-Weil bound**](https://www.jstor.org/stable/pdf/2372655.pdf). Giving all the details goes to far, but the cartoon-version is: <br>\n",
    "The space $\\mathbb{F}_q^3$ is 3-dimensional, and the solution set $\\mathbb{L}_q(f) \\subseteq \\mathbb{F}_q^3$ is cut out by $1$ equation, and has dimension $3-1=2$. It is also *irreducible*, i.e. it is not (for all $q$) the union of the solution sets of two equations $f_1, f_2$. Then the Lang-Weil bound says:\n",
    "$$\n",
    "|\\mathbb{L}_q(f)| = q^{\\dim \\mathbb{L}_q(f)} + O(q^{\\dim \\mathbb{L}_q(f)-1/2}) = q^2 + O(q^{1.5})\n",
    "$$\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Symbolic expressions\n",
    "We have already seen in some examples that SageMath can do computations involving formal variables. These happen in the so-called **Symbolic Ring**, which is denoted by ``SR`` in SageMath. \n",
    "\n",
    "### Creating and comparing symbolic expressions\n",
    "To get started, we declare some variable names using the function ``var``."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('a,b')\n",
    "f = (a+b)^3; f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(f)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f in SR"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Then we can ask to perform operations, like expanding the formula we entered above:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "expand(f)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also try to check whether two formulas are equivalent:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "t1 = (a^2+2*a*b+b^2)/(a+b); t1"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "t2 = a+b; t2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "t1 == t2"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bool(t1 == t2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "This is not restricted to algebraic equations:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "bool(sin(a)^2 + cos(a)^2 == 1)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Try to check the following identities, some of which are true and some of which are false. Which one does SageMath get wrong?\n",
    "$$\n",
    "a \\leq \\max(a,b)\n",
    "\\\\\n",
    "\\exp(a) \\cdot \\exp(b) = \\exp(a+b)\n",
    "\\\\\n",
    "\\frac{1}{\\frac{1}{a}-\\frac{1}{b}} = \\frac{a \\cdot b}{a - b}\n",
    "\\\\\n",
    "\\frac{1}{2}(a+b) \\leq \\max(a,b)\n",
    "\\\\\n",
    "(a+b+c)^2 = a^2 + b^2 + c^2 + 2 \\cdot(ab + ac + bc)\n",
    "\\\\\n",
    "\\sin\\left(\\frac{\\pi}{4}\\right) = \\frac{1}{\\sqrt{2}}\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)<br>\n",
    "<!---\n",
    "The first two statements are correct, and are evaluated as True:\n",
    "```\n",
    "bool(a <= max(a,b))\n",
    "> True\n",
    "bool(exp(a)*exp(b)==exp(a+b))\n",
    "> True\n",
    "```\n",
    "The next statement is false (there is a sign mistake in the denominator):\n",
    "```\n",
    "bool(1/(1/a - 1/b) == a*b/(a-b))\n",
    "> False\n",
    "```\n",
    "The inequality\n",
    "$$\n",
    "\\frac{1}{2}(a+b) \\leq \\max(a,b)\n",
    "$$\n",
    "is correct, but unfortunately SageMath evaluates it as ``False``:\n",
    "```\n",
    "bool(1/2*(a+b) <= max(a,b))\n",
    "> False\n",
    "```\n",
    "One should note that [other computer algebra systems get this right](https://www.wolframalpha.com/input?i=1%2F2*%28a%2Bb%29%3C%3Dmax%28a%2Cb%29).\n",
    "\n",
    "The last two formulas are again correct. Since we need an additional variable, we create it using ``var``:\n",
    "```\n",
    "var('c')\n",
    "bool((a+b+c)^2 == a^2 + b^2 + c^2 + 2*(a*b + a*c + b*c))\n",
    "> True\n",
    "bool(sin(pi/4) == 1/sqrt(2))\n",
    "> True\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Substitution\n",
    "We can make variable substitutions as follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('a,b,c,d')\n",
    "f = a^2 + b*c; f"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f.subs({a:b, c:b+1})"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Note that the argument of ``f.subs`` is a dictionary, sending variables appearing in ``f`` to their new values. Alternatively, we can treat ``f`` as if it was a function, and call it (but specifying the names of the variables as follows):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(a=2,b=3,c=4)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "f(a=d,b=d,c=d)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Formal sums\n",
    "We can also evaluate some formal summations, e.g. the famous identity\n",
    "$$\n",
    "1 +2 + \\ldots + n = \\frac{n(n+1)}{2}\n",
    "$$\n",
    "Indeed, to compute\n",
    "$$\n",
    "\\sum_{a=a_0}^{a_1} f(a)\n",
    "$$\n",
    "we use ``sum(f,a,a0,a1)``, where ``a`` is a variable and ``f`` is again a symbolic expression."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('n')\n",
    "g = sum(a,a,1,n); g"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The answer is again a symbolic expression, so we can e.g. substitute some value for ``n`` or compute a further symbolic sum:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g.subs({n:100}) # The famous sum that Gauss allegedly computed in elementary school"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('m')\n",
    "sum(g,n,1,m)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### Solving equations\n",
    "Maybe you have wondered why writing ``f == g`` for two symbolic expressions does not immediately give either ``True`` or ``False``, but again a symbolic expression:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "g == 5050"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "type(g == 5050)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The reason is that we can now use this expression to state an equation (or a system of equations) that we want to solve:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solve(g == 5050, n)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The second argument is the variable for which we want to solve. To specify multiple equations (or inequalities) we give a list of these as the first argument of ``solve``:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solve([g == 5050, n>=0], n)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To get an output that can be used in further computations, we can specify the optional parameter ``solution_dict=True`` to obtain a list of dictionaries giving our solutions:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solve([g == 5050, n>=0], n, solution_dict=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('a,b')\n",
    "L = solve(a^2 - b^2 == 0, b, solution_dict=True); L"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "for l in L:\n",
    "    print((a^2 - b^2).subs(l))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Finally, we can also solve for multiple variables:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "var('y')\n",
    "solve([x*y==0, (x-1)*(y-1)==0],[x,y])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Solve the following (systems of) equalities for the variables $x,y,n,\\ldots$:\n",
    "$$\n",
    "x^3 + 2 x^2 - 21 x + 18 = 0, \\\\\n",
    "x^2 + p x + q = 0,\\\\\n",
    "\\sum_{k=0}^n 2^k = 4095\\\\\n",
    "x^2 + y^2 = 1, x + y = 1\n",
    "$$"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)<br>\n",
    "<!---\n",
    "Before we start, we should introduce all variables that we need:\n",
    "```\n",
    "var('x,y,p,q,k,n')\n",
    "```\n",
    "$$\n",
    "x^3 + 2 x^2 - 21 x + 18 = 0\n",
    "$$\n",
    "```\n",
    "solve(x^3 + 2*x^2 - 21*x + 18 == 0, x)\n",
    "> [n == 11]\n",
    "```\n",
    "$$\n",
    "x^2 + p x + q = 0\n",
    "$$\n",
    "```\n",
    "solve(x^2 + p*x + q == 0, x)\n",
    "> [x == -1/2*p - 1/2*sqrt(p^2 - 4*q), x == -1/2*p + 1/2*sqrt(p^2 - 4*q)]\n",
    "```\n",
    "$$\n",
    "\\sum_{k=0}^n 2^k = 4095\n",
    "$$\n",
    "```\n",
    "solve(sum(2^k,k,0,n)==4095,n)\n",
    "> [n == 11]\n",
    "```\n",
    "$$\n",
    "x^2 + y^2 = 1, x + y = 1\n",
    "$$\n",
    "```\n",
    "solve([x^2 + y^2 == 1, x+y == 1], [x,y])\n",
    "> [[x == 1, y == 0], [x == 0, y == 1]]\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We can also add some [assumptions](https://doc.sagemath.org/html/en/reference/calculus/sage/symbolic/assumptions.html) about the variables:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "solve(x^3==1, x, solution_dict=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "assume(x, 'real')\n",
    "solve(x^3==1, x, solution_dict=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Here is how we get rid of the assumptions again:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "forget(assumptions())\n",
    "solve(x^3==1, x, solution_dict=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Assignments"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Given a number $d \\in \\{2,3, \\ldots, 9\\}$ here are the rules of a game named $d$-Blup:\n",
    "\n",
    "The participants stand in a circle and count upwards starting from $1$, except that for every number either divisible by $d$ or ending on the digit $d$, they have to say \"Blup\" instead. If you say something wrong or take too long, you are out, and the others start again. \n",
    "\n",
    "A correct $3$-Blup would start as follows:\n",
    "```\n",
    "Player  1 : \"1\"\n",
    "Player  2 : \"2\"\n",
    "Player  3 : \"Blup\"\n",
    "Player  4 : \"4\"\n",
    "...\n",
    "Player 11 : \"11\"\n",
    "Player 12 : \"Blup\"\n",
    "Player 13 : \"Blup\"\n",
    "...\n",
    "```\n",
    "Given ``d``, write a program that prints the first 50 things the players need to say, i.e.\n",
    "```\n",
    "1\n",
    "2\n",
    "Blup\n",
    "4\n",
    "...\n",
    "```\n",
    "*Hint:* Recall that given integers ``a,b`` with ``b`` nonzero, you can compute the remainder of the division of ``a`` by ``b`` using ``a % b``."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution:** (uncomment to see)<br>\n",
    "<!---\n",
    "Below you see code which does the trick (note that ``n % 10`` gives the last digit of ``n``). The line ``sleep(1)`` makes the program wait ``1`` second, so you can play the game against the computer.\n",
    "```\n",
    "for n in range(1,51):\n",
    "    sleep(1)\n",
    "    if n % d == 0 or n % 10 == d:\n",
    "        print('Blup')\n",
    "    else:\n",
    "        print(n)\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "The [Collatz conjecture](https://en.wikipedia.org/wiki/Collatz_conjecture) states that for the function\n",
    "$$\n",
    "f : \\mathbb{Z}_{\\geq 1} \\to \\mathbb{Z}_{\\geq 1}, n \\mapsto \n",
    "\\begin{cases}\n",
    "3n+1 & \\text{for $n$ odd},\\\\\n",
    "n/2 & \\text{for $n$ even}\n",
    "\\end{cases}\n",
    "$$\n",
    "and any positive integer $n$, the sequence $n$, $f(n)$, $f(f(n))$, $\\ldots$ eventually hits the number $1$. Test the Collatz conjecture for the first $100$ integers.<br>\n",
    "*Hint:* The ``while``-loop\n",
    "```\n",
    "while CONDITION:\n",
    "    do_something()\n",
    "```\n",
    "repeats the following code block until the ``CONDITION`` becomes ``False``."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)\n",
    "<!---\n",
    "```\n",
    "for n in range(1,101):\n",
    "    N = n\n",
    "    while(N != 1):\n",
    "        if N % 2 == 0:\n",
    "            N = N/2\n",
    "        else:\n",
    "            N = 3*N+1\n",
    "    print('{} is no counterexample to the Collatz conjecture.'.format(n))\n",
    "```\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "The Pólya conjecture states that given a natural number $n>1$, there are among the set\n",
    "$$\\{1, \\ldots, n\\}$$\n",
    "at least as many elements with an odd number of prime factors as elements with an even number of prime factors. Here prime factors are counted with multiplicity, i.e. $12 = 2 \\cdot 2 \\cdot 3$ has $3$ prime factors according to this definition.\n",
    "* Test the Pólya conjecture for the integers $2 \\leq n \\leq n_0$, for an $n_0$ of your choice.\n",
    "* Google the Pólya conjecture."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution:** (uncomment to see)<br>\n",
    "<!---\n",
    "Below we go through the numbers $n=2,3, \\ldots$ and have variables ``odd_pf`` and ``even_pf`` which record how many numbers we have already seen that have odd and even numbers of prime factors. If ever the number ``even_pf`` becomes strictly bigger than ``odd_pf``, we found a counter-example and print it, breaking the for loop with the command ``break`` (this is optional).\n",
    "```\n",
    "odd_pf = 0\n",
    "even_pf = 1\n",
    "for n in range(2,10000000):\n",
    "    factors = dict(factor(n))\n",
    "    if sum(factors.values())%2 == 1:\n",
    "        odd_pf += 1\n",
    "    else:\n",
    "        even_pf += 1\n",
    "    if even_pf > odd_pf:\n",
    "        print('Counterexample:')\n",
    "        print(n)\n",
    "        break\n",
    "```\n",
    "After consulting the [wikipedia page of the conjecture](https://en.wikipedia.org/wiki/P%C3%B3lya_conjecture) you see that you need some patience to find the first counter-example $n = 906,150,257$, which took my computer about 4 hours.\n",
    "\n",
    "**Challenge exercise**<br>\n",
    "Write a program which finds this counter-example in less than one hour.<br>\n",
    "*Disclaimer:* I have not yet solved this exercise myself yet, which is why it's a *challenge* exercise ...\n",
    "--->"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "#### Exercise\n",
    "Given an integer $c \\geq 0$ let $T_c \\subseteq \\mathbb{Z}^2$ be the set of points in $\\mathbb{Z}^2$ contained in the triangle spanned by $(0,0), (0,c), (c,0)$. Or in other words\n",
    "$$\n",
    "T_c = \\{(x,y) \\in \\mathbb{Z}^2 : x,y \\geq 0, x+y \\leq c\\}.\n",
    "$$\n",
    "We claim that the function\n",
    "$$\n",
    "g : \\mathbb{Z}_{\\geq 0} \\to \\mathbb{Z}, c \\mapsto \\sum_{(x,y) \\in T_c} x \\cdot y\n",
    "$$\n",
    "is given by a polynomial in $c$. Compute this polynomial and check your result in a few cases."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "**Solution** (uncomment to see)<br>\n",
    "<!---\n",
    "We can write the function $g$ as\n",
    "$$\n",
    "g(c) = \\sum_{x=0}^c \\sum_{y=0}^{c-x} x \\cdot y.\n",
    "$$\n",
    "By setting\n",
    "$$\n",
    "f(x,c) = \\sum_{y=0}^{c-x} x \\cdot y\n",
    "$$\n",
    "it is now very easy to compute this using symbolic sums:\n",
    "```\n",
    "var('x,y,c')\n",
    "f = sum(x*y,y,0,c-x)\n",
    "g = sum(f,x,0,c); g\n",
    "> 1/24*c^4 + 1/12*c^3 - 1/24*c^2 - 1/12*c\n",
    "```\n",
    "We can now compute a value, e.g. at $c=2$ where we have\n",
    "$$\n",
    "g(2) = 0 \\cdot 0 + 1 \\cdot 0 + 0 \\cdot 1 + 2 \\cdot 0 + 1 \\cdot 1 + 0 \\cdot 2 = 1,\n",
    "$$\n",
    "which we verify as follows:\n",
    "```\n",
    "g.subs({c:2})\n",
    "> 1\n",
    "g.subs({c:3})\n",
    "> 5\n",
    "```\n",
    "The value at $c=3$ is explained by the additional summands $2 \\cdot 1 + 1 \\cdot 2 = 4$.\n",
    "\n",
    "To prove the formula you can use [Faulhaber's formulas](https://en.wikipedia.org/wiki/Faulhaber%27s_formula#Examples) for the sum of the $k$-th powers of the first $n$ integers. More generally, the fact that this summation over a triangle is a polynomial follows from a variant of [Ehrhart theory](https://en.wikipedia.org/wiki/Ehrhart_polynomial).\n",
    "--->"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "SageMath 9.1",
   "language": "sage",
   "name": "sagemath"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.7.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}